Probabilistic Seismic Hazard Analysis Research Articles

A prolegomenon to Design Input Energy Spectra for the Himalayan region

Energy-based design is becoming progressively significant as a new age of earthquake-resistant design unfolds. Thitherto, very little is known about energy-based design in India despite the fact that the Himalayan region is one of the most seismically active regions in the world. Therefore, the current work provides a preliminary analysis of the elastic input energy and equivalent velocity (Vea) spectra for the active regions in the Himalayas. In this regard, the ground motions of Himalayan databases are blended with NGA WEST2 to estimate the Vea spectra, an inertia-independent representation of the input energy spectra. Incipiently, an Artificial Neural Network-based ground motion model (GMM) is developed for Vea spectra based on the predictor variables as magnitude, distance, hypocentral depth and site conditions. Consequently, a Probabilistic Seismic Hazard analysis is performed to obtain the uniform hazard energy spectra for a 2475-year return period. In addition, a Design Input Energy Spectra is proposed for the active regions in the Indian Himalayas. In addition, the relationship between the Vea and pseudo-spectral velocity (PSV) is analyzed to develop a Vea-informed GMM for PSV. Furthermore, the PSA estimated from the Vea-PSV model is compared with the design spectra recommended by the Indian Standard Code of Practice (IS1893-1:2016).

Seismic Hazard Microzonation Map for the Central Plain of Thailand

A study was conducted to investigate the seismic hazard microzonation of the central plain of Thailand, which is situated in an area with a thick quaternary basin. Unconsolidated sediments can lead to amplification of earthquake ground shaking at fundamental frequency and can cause a significant increase in damage to buildings. The research yielded significant results, including the development of a fundamental frequency map through the analysis of HVSR for 149 microtremor measurement sites. Subsequently, a Vs30 map was derived utilizing HVSR inversion techniques, and a soil classification map was constructed based on the NEHRP classification. The upper central plain along the Yom River and the Nan River had a low fundamental frequency of 0.3-0.5 Hz and low Vs30, which can be classified as soil type Class E. In the southern areas of Ayutthaya, Pathum Thani, and central Bangkok, an extremely low Vs30 of less than 100 m/s was observed, indicating soil class F or special soft soil. A comprehensive investigation was conducted on probabilistic seismic hazard analysis by considering the Vs30 sites condition for PGA, SA0.2s, and SA1.0s with a 2475-year return period. The northern region of the upper central plain and the western side of the central basin exhibit relatively high seismic hazards. Furthermore, the site effect significantly amplifies ground motion at the 1.0 second period, surpassing the earthquake-resistant design standard for buildings in Thailand by more than 5 times, particularly in the central region of the lower central plain.

Regional and Site-Specific Ground Motion Model for Probabilistic Seismic Hazard Analysis in Taiwan: A Case Study of I-Lan

In this study, a probabilistic seismic hazard analysis (PSHA) using a traditional approach with a simple logic tree that handles a single-site sigma is developed to present a site-specific ground motion hazard. The analysis relies on the accelerometer data collected from the I-Lan Plain, a deep sedimentary basin in Northeastern Taiwan where hundreds of earthquake sequences (with [Formula: see text] = 4.0–7.6) were recorded by the Taiwan Strong Motion Instrument Program (TSMIP) network between 1992 and 2016. The performance of PSHA on soil and hard sites is evaluated from two representative instrumented sites (ILA048 and ILA025), using residual measures available from a ground motion model (GMM) developed by Phung et al. [ 2023 ] (referred to as Ph23), from which the uncertainties in site terms ([Formula: see text] and single-site sigma ([Formula: see text] are estimated. We address key conceptual issues in the current PSHA approach and introduce a new region- and site-specific PSHA approach in which (1) site-to-site variability ([Formula: see text] is estimated as a random variance in a mixed effects GMM regression and (2) the GMM site-specific single-site sigma ([Formula: see text] is replaced with a generic site-corrected aleatory variability ([Formula: see text]). Comparison of the region- and site-specific hazard curves from our method against the traditional method estimates at two well-recorded sites in the I-Lan region shows an approximate 50% difference in prediction ground motion values considering for 2% and 10% probability of exceedance in 50 years.

Evaluation of a Decade-Long Prospective Earthquake Forecasting Experiment in Italy

Abstract Earthquake forecasting models represent our current understanding of the physics and statistics that govern earthquake occurrence processes. Providing such forecasts as falsifiable statements can help us assess a model’s hypothesis to be, at the least, a plausible conjecture to explain the observations. Prospective testing (i.e., with future data, once the model and experiment have been fully specified) is fundamental in science because it enables confronting a model with completely out-of-sample data and zero degrees of freedom. Testing can also help inform decisions regarding the selection of models, data types, or procedures in practical applications, such as Probabilistic Seismic Hazard Analysis. In 2010, a 10-year earthquake forecasting experiment began in Italy, where researchers collectively agreed on authoritative data sources, testing rules, and formats to independently evaluate a collection of forecasting models. Here, we test these models with ten years of fully prospective data using a multiscore approach to (1) identify the model features that correlate with data-consistent or -inconsistent forecasts; (2) evaluate the stability of the experiment results over time; and (3) quantify the models’ limitations to generate spatial forecasts consistent with earthquake clustering. As each testing metric analyzes only limited properties of a forecast, the proposed synoptic analysis using multiple scores allows drawing more robust conclusions. Our results show that the best-performing models use catalogs that span over 100 yr and incorporate fault information, demonstrating and quantifying the value of these data types. Model rankings are stable over time, suggesting that a 10-year period in Italy can provide sufficient data to discriminate between optimal and suboptimal forecasts. Finally, no model can adequately describe spatial clustering, but those including fault information are less inconsistent with the observations. Prospective testing assesses relevant assumptions and hypotheses of earthquake processes truly out-of-sample, thus guiding model development and decision-making to improve society’s earthquake resilience.

The Observed and Simulated Site Effects at the Selected Regions Influenced by 2023 Turkey Earthquake Sequences: A Preliminary Evaluation

ABSTRACT This study offers an initial evaluation of site effects within regions impacted by the February 2023 earthquakes in Turkey. It develops a preliminary response spectrum site amplification model, alongside probabilistic seismic hazard analysis and 1D site response analyses (SRAs) specifically focused on a site in Hatay. The model suggests higher site factors in softer sites for 0.2 sec, contrasting with other conditions. Results demonstrate that surface spectra derived from 1D SRAs and site-specific uniform hazard spectrum are consistently lower than Turkish Earthquake Building Code (TEBC) design spectrum across all considered periods, indicating the necessity for a detailed reassessment of TEBC specifications.

The Upper Crustal Deformation Field of Greece Inferred From GPS Data and Its Correlation With Earthquake Occurrence

AbstractWe present a new geodetic strain rate and rotation rate model for Greece that has been derived using a highly dense GPS velocity field. The spatial distribution and the resolved rates of the various velocity gradient tensor quantities provided updated constraints on the present‐day upper crustal deformation in the region and revealed new details not reported previously. The spatial distribution of the second invariant demonstrated that the overall magnitude of strain rates is highest across two well‐defined provinces. The first follows the North Anatolian Fault and its two branches within the north Aegean, crosses central Greece and through the Gulf of Corinth it terminates in western Greece, while the second encompasses the extensional province of western Turkey and the eastern Aegean Sea islands. Our estimates revealed that shearing affects some of the fault‐bounded grabens of central Greece that lie to the SW of the North Aegean Basin implying considerable oblique extension. We identified a narrow region of counterclockwise rotation whose location and kinematics have been induced by the net effect across the intersection of the clockwise rotating domains of western and central Greece. The Aegean microplate and the Anatolian plate are separated by a wide transition zone which accommodates the curved stretching of the entire plate system. In both edges of the Hellenic forearc the dominant mode of crustal strain is E‐W extension. We found that earthquakes of M ≥ 5.6 are spatially well‐correlated with high‐strain areas, indicating that strain rate mapping could be used to inform future probabilistic seismic hazard analyses.

Prediction of Fourier Amplitude Spectrum of Ground Motion in Mexico City from Subduction Thrust Earthquakes

Mexico City, one of the largest cities in the world, experiences frequent, devastating interplate earthquakes that originate on the subduction thrust along the Pacific coast of Mexico more than 250 km away. A notorious example is the 19 September 1985 Michoacán earthquake (Mw 8.0) which killed ~ 10,000 persons and caused wide-spread destruction in the city. The main cause of seismic damage in Mexico City is its subsoil that amplifies the ground motion. In view of the seismic hazard faced by the city, a reliable estimation of ground motion during future earthquakes is of vital importance. We present a new ground motion prediction equation (GMPE) for the Fourier amplitude acceleration spectrum (FAS) from subduction thrust earthquakes along the Pacific coast of Mexico at the reference station, CU, located in the firm zone of Mexico City. The GMPE is derived via Bayesian regression analysis and is based on an enlarged set of recordings (1965-2020; 40 earthquakes; 250 ≤ R≤ 500 km; 5 ≤ Mw ≤ 8.0). An important feature of the new GMPE is that it includes a term to model different attenuation along different ray paths. The inclusion of this term leads to a reduction in the aleatory variability of the GMPE, particularly at high frequencies. Since spectral amplification of seismic waves with respect to CU is known at many sites in the city, the FAS at these sites can be computed if it is known at CU. The new GMPE has been incorporated in a fully Fourier-based probabilistic seismic hazard analysis of Mexico City.

The influence of aftershocks on seismic hazard analysis: a case study from Xichang and the surrounding areas

Abstract. In some instances, a strong aftershock can cause more damage than the mainshock. Ignoring the influence of aftershocks may lead to the underestimation of the seismic hazard of some areas. Taking Xichang and its surrounding areas as an example and based on the Seismic ground motion parameters zonation map of China (GB 18306-2015), this study used the Monte Carlo method to simulate synthetic mainshock sequences. Additionally, the Omi–Reasenberg–Jones (Omi–R–J) aftershock activity model is used to simulate the aftershock sequences that follow mainshocks above a certain magnitude threshold. Then, the mainshock and the aftershocks are combined to calculate the regional seismic hazard using ground motion prediction equations (GMPEs). Finally, the influence of aftershocks on seismic hazard analysis is examined and considered. The results show that in areas with moderate to strong seismic backgrounds, the influence of aftershocks on probabilistic seismic hazard analysis can exceed 50 %. These results suggest that the impact of aftershocks should be properly considered for future probabilistic seismic hazard analyses, especially in areas with moderate to strong seismic activity backgrounds and in areas prone to secondary disasters such as landslides and mudslides.

Integration of rupture directivity models for the US National Seismic Hazard Model

Several rupture directivity models (DMs) have been developed in recent years to describe the near-source spatial variations in ground-motion amplitudes related to propagation of rupture along the fault. We recently organized an effort toward incorporating these directivity effects into the US Geological Survey (USGS) National Seismic Hazard Model (NSHM), by first evaluating the community’s work and potential methods to implement directivity adjustments into probabilistic seismic hazard analysis (PSHA). Guided by this evaluation and comparison among the considered DMs, we selected an approach that can be readily implemented into the USGS hazard software, which provides an azimuthally varying adjustment to the median ground motion and its aleatory variability. This method allows assessment of the impact on hazard levels and provides a platform to test the DM amplification predictions using a generalized coordinate system, necessary for consistent calculation of source-to-site distance terms for complex ruptures. We give examples of the directivity-related impact on hazard, progressing from a simple, hypothetical rupture, to more complex fault systems, composed of multiple rupture segments and sources. The directivity adjustments were constrained to strike–slip faulting, where DMs have good agreement. We find that rupture directivity adjustments using a simple median and aleatory adjustment approach can affect hazard both from a site-specific perspective and on a regional scale, increasing ground motions off the end of the fault trace up to 30%–40% and potentially reducing it for sites along strike. Statewide hazard maps of California show that the change in shaking along major faults can be a factor to consider for assessing long-period [Formula: see text] near-source effects within the USGS NSHM going forward, reaching up to 10%–20%. Finally, we suggest consideration of minimum parameter ranges and baseline requirements as future DMs are developed to minimize single approach adaptations to enable more consistent application within both ground motion and hazard studies.

Site-specific stochastic ground motion model utilizing deterministic physics-informed simulations: A Bayesian approach

Limited availability of recorded ground motions poses a challenge for reliable probabilistic seismic-hazard analysis (PSHA), even in highly seismic regions like the Western United States. Stochastic ground motions are commonly employed to address this challenge. However, the stochastic ground motion models (GMMs) may not consistently generate ground motions compatible with the site hazard due to their calibration using global data, failing to capture site-specific characteristics adequately. In the absence of recorded motions, physics-informed simulations provide a viable alternative but are deterministic with limitations of their own that makes them challenging to support PSHA. This article introduces a Bayesian framework that combines prior knowledge from a stochastic GMM, calibrated with global data, with site-specific data obtained from deterministic physics-informed simulations. The proposed framework utilizes the Rezaeian–Der Kiureghian (2010) model as the stochastic GMM and incorporates site-specific data from the CyberShake 15.12 study. By updating the mean and variance of the predictive relationships, along with the marginal distribution of the model parameters, through Bayesian inference, this framework allows for the simulation of site-specific ground motions consistent with the site characteristics. The statistics of peak ground acceleration distributions, as well as both the median and variability of the elastic response spectra, obtained from the calibrated stochastic GMM, demonstrate consistency with those derived using GMMs based on the Next Generation Attenuation (NGA) database.

Study on seismic risk assessment model of water supply systems in mainland China

Abstract. Using the peak ground acceleration (PGA) under four exceedance probabilities calculated via a probabilistic seismic hazard analysis method for China, the probability density function of PGA was obtained by fitting with the Cornell seismic hazard exceedance probability PGA function model. Combined with the seismic fragility function of the water supply system calculated on the basis of the empirical matrix of actual earthquake damage and the exposure of fixed assets, the expected seismic loss and expected loss rate models of the water supply system were established, and the classification standard with the expected seismic loss rate of the water supply system as the index was proposed. The seismic fragility of the water supply system was classified, and the exposure of the water supply system was analyzed. The total fixed assets in the Statistical Yearbook of Urban Water Supply, henceforth called Water Supply Yearbook, were taken as the exposure to earthquake in the region. The accuracy of the fragility model in this paper was verified through the actual earthquake damage losses in Deyang. Taking the water supply system of 720 cities in mainland China as an example, the distribution maps of expected seismic loss and expected loss rate were calculated and drawn. The expected loss rate model was verified by the key earthquake prevention areas in mainland China. The assessment model based on the expected loss and expected loss rate was taken as the seismic risk assessment model of water supply systems in mainland China.

The 2023 US National Seismic Hazard Model: Ground-motion characterization for the conterminous United States

We update the ground-motion characterization for the 2023 National Seismic Hazard Model (NSHM) for the conterminous United States. The update includes the use of new ground-motion models (GMMs) in the Cascadia subduction zone; an adjustment to the central and eastern United States (CEUS) GMMs to reduce misfits with observed data; an updated boundary for the application of GMMs for shallow, crustal earthquakes in active tectonic regions (i.e. western United States (WUS)) and stable continental regions (i.e. CEUS); and the use of improved models for the site response of deep sedimentary basins in the WUS and CEUS. Site response updates include basin models for the California Great Valley and for the Portland and Tualatin basins, Oregon, as well as long-period basin effects from three-dimensional simulations in the Greater Los Angeles region and in the Seattle basin; in the CEUS, we introduce a broadband (0.01- to 10-s period) amplification model for the effects of the passive-margin basins of the Atlantic and Gulf Coastal Plains. In addition, we summarize progress on implementing rupture directivity models into seismic hazard models, although they are not incorporated in the 2023 NSHM. We implement the ground-motion characterization for the 2023 NSHM in the US Geological Survey’s code for probabilistic seismic hazard analysis, nshmp-haz-v2, and present the sensitivity of hazard to these changes. Hazard calculations indicate widespread effects from adjustments to the CEUS GMMs, from the incorporation of Coastal Plain amplification effects, and from the treatment of shallow-basin and out-of-basin sites in the San Francisco Bay Area and Los Angeles region, as well as locally important changes from subduction-zone GMMs, and from updated and new WUS basins.

Territorial exceedance of probabilistic seismic hazard from ShakeMap data

Current seismic structural design makes use of a ground motion intensity that has a certain probability of being exceeded at a site of interest in a time interval or, equivalently, exceedance return period. The design intensities with the same return period are often collected in the form of maps deriving from probabilistic seismic hazard analysis (PSHA) for each of the sites of interest. Probability theory underlying PSHA dictates that, in any time interval, design intensities are expected to be exceeded in a fraction of sites that depends on the return period the map refers to. In the case of Italy, three different nationwide PSHA studies can be currently considered of relevance. In the study, the estimated areal fraction of the Italian territory exposed to exceedance of the design intensity from 2008 to 2019 was quantified for the three hazard models, based on ShakeMap data for instrumental earthquakes. In addition, the same fraction was calculated considering a sparse catalog of inferred ShakeMap for historical earthquakes that occurred over almost 1000 years. It was found that, despite the apparent differences in the hazard models, the estimated fraction of territory exposed to exceedance is comparable for all the considered hazard maps.

Analysis of Soil Vulnerability Level Due to Earthquake in Tarutung Region and Its Surroundings Using Probabilistic Seismic Hazard Analysis Method

Analysis of Soil Vulnerability Level Due to Earthquake in Tarutung Region and Its Surroundings Using Probabilistic Seismic Hazard Analysis Method

Uncertainty and sensitivity assessments on the inputs of probabilistic seismic hazard assessment: A case study of the North Tehran fault

Probabilistic seismic hazard analysis (PSHA) often incorporates significant uncertainty in its inputs, which have the potential to substantially affect the estimated hazards. To improve the reliability of seismic hazards, the present work examines how these uncertainties impact the variability of PSHA and disaggregation results. For this purpose, uncertainty and sensitivity analyses have been performed for different sites in the Tehran region, with the seismic source of the North Tehran Fault as a case study. The results show that uncertainties in the input parameters have a moderate to significant effect on the PSHA and disaggregation results. The uncertainty and sensitivity analyses of this study provide a better understanding of the various PSHA parameters and their associated uncertainties. Therefore, such analyses establish confidence in the PSHA results and expose the main sources of uncertainty, which enable the prioritization of endeavors to reduce hazard variability.

PSHA: Does It Deal with What It Is or What We Want It to Be?

Abstract In a recent opinion piece Albarello and Paolucci (2023; hereafter, AP23) provide their view as members of the past Seismic Group of the Commissione Grandi Rischi (CGR-SRS) in Italy, which represents the main scientific consultant for Italian Civil Protection, about the difficulty using probabilistic seismic hazard analysis (PSHA) models for building code purposes. Here, we refer to this specific kind of PSHA modeling as National Seismic Hazard Model (NSHM). We agree with AP23 that the topic is of great and general importance, and here we aim at contributing to this discussion by offering our perspective on two points that are at the heart of the matter, concluding that AP23 is misguided in how to deal with them. First, we assert that the credibility of an NSHM has to be rooted only in the use of the best available science, which includes a rigorous testing phase with observations, independent from the consequences in terms of risk. (PSHA deals with what it is.) Second, we claim that the difficulties in accepting a new NSHM with some major changes with respect to the previous model are mostly due to too rigid building code procedures that do not account for the epistemic uncertainty in the hazard estimates.

Aftershock probabilistic seismic hazard analysis based on enhanced Bayesian network considering frequency information

AbstractBayesian network (BN) is an important tool in probabilistic seismic risk analysis (PSRA) due to its holistic nature and powerful probabilistic inference capabilities. However, while the information that can be stored by BN includes variable values, probability distributions, and relationships between variables, it does not involve event quantities. The results obtained from PSRA and probabilistic seismic hazard analysis (PSHA) are commonly frequency‐related (the number of occurrences per unit time), fundamentally describing event numbers. BN assumes a fixed event number relationship between nodes, but the aftershock counts vary for different mainshock magnitudes. This limitation further restricts its effectiveness in lifespan PSRA for structures. Therefore, this study aims to propose an enhanced BN (EBN) for modeling aftershock PSHA (APSHA) considering frequency information. This paper first introduces the modeling method using BN, including PSHA and a simplified model of APSHA using Båth's law. Next, an EBN that can effectively address frequency issues between the mainshock and aftershocks is presented for modeling APSHA considering the modified Omori's law, while the mathematical principles of both forward and backward inferences for the EBN are elaborated in detail. Then, the article carefully analyzes the utilization of the EBN for APSHA modeling under different scenarios where the information regarding the mainshock is known to varying extents, followed by the extension of the EBN in a time domain to obtain time‐dependent aftershock hazard. Finally, an example of APSHA in a specific location in China is illustrated.

Comparison of peak ground acceleration using deterministic and probabilistic approach in the coastal area of Loh Buaya, Rinca Island, Indonesia

An area located near an active volcano, subduction zone, and/or active fault has a risk of an earthquake. Earthquakes can cause damage to buildings with a large number of losses. One of the earthquake-related parameters required in planning an earthquake-resistant building is the Peak Ground Acceleration (PGA). The objective of this study is to compare the PGA value at the coastal area of Loh Buaya, Rinca Island, Indonesia, using Probabilistic Seismic Hazard Analysis (PSHA) and Deterministic Seismic Hazard Analysis (DSHA). The DSHA is calculated based on earthquake data series from the United States Geological Survey (USGS) in the period between 1922 to 2022. Also, it is calculated using 10 data of Ground Motion Prediction Equation (GMPE). Data from the USGS is constrained with magnitude > 5.0 MW and earthquake distance < 200 km from the study site. The PSHA is based on the Indonesian earthquake hazard map and the de-aggregation map from the National Earthquake Study Centre of Indonesia with a probability exceeding 2% in 50 years. Based on the parameters of PGA maximum, distance, and magnitude, the DSHA by Kanno method has the closest result to the PSHA. It is concluded that this area has a high earthquake risk with a PGA maximum of more than 0.5 g.

MCPSHA: A New Tool for Probabilistic Seismic Hazard Analysis Based on Monte Carlo Simulation

The utilization of the Monte Carlo method in conjunction with probabilistic seismic hazard analysis (PSHA) constitutes a compelling avenue for exploration. This approach presents itself as an efficient and adaptable alternative to conventional PSHA, particularly when confronted with intricate factors such as parameter uncertainties and diverse earthquake source models. Leveraging the Monte Carlo method and drawing from the widely adopted Cornell-type seismicity model in engineering seismology and disaster mitigation, as well as a seismicity model capturing temporal, spatial, and magnitude inhomogeneity, we have derived a formula for the probability of earthquake intensity occurrence and the mean rate of intensity occurrence over a specified time period. This effort has culminated in the development of a MATLAB-based program named MCPSHA. To assess the model's efficacy, we selected Baoji City, Shaanxi Province, China, as our research site. Our investigation delves into the disparity between occurrence probability and extreme probability (a surrogate commonly employed for occurrence probability) in the Baoji region over the next 50 years. The findings reveal that the Western region of Baoji exhibits a heightened hazard level, as depicted in the maps, which illustrate a 10% probability of exceedance within a 50-year timeframe. The probability of earthquake occurrence under various intensities (VI, VII, and VIII) over 50 years follows a declining trend from west to east. Furthermore, the likelihood of seismic intensity exceeding VI, VII, and VIII indicates the lowest exceeding probability in the northeast and the highest in the northwest. Notably, for intensities VI-VII, the difference between occurrence probability and extreme probability approaches twice, gradually diminishing with increasing intensity. This study underscores the MCPSHA model's efficacy in providing robust technical support for mitigating earthquake risk and enhancing the precision of earthquake insurance premium rate calculations.

Accounting for the Variability of Earthquake Rates within Low-Seismicity Regions: Application to the 2022 Aotearoa New Zealand National Seismic Hazard Model

ABSTRACT The distribution of earthquakes in time and space is seldom stationary, which could hinder a robust statistical analysis, particularly in low-seismicity regions with limited data. This work investigates the performance of stationary Poisson and spatially precise forecasts, such as smoothed seismicity models (SSMs), in terms of the available training data. Catalog bootstrap experiments are conducted to: (1) identify the number of training data necessary for SSMs to perform spatially better than the least-informative Uniform Rate Zone (URZ) models; and (2) describe the rate temporal variability accounting for the overdispersion and nonstationarity of seismicity. Formally, the strict-stationarity assumption used in traditional forecasts is relaxed into local and incremental stationarity (i.e., a catalog is only stationary in the vicinity of a given time point t) along with self-similar behavior described by a power law. The results reveal rate dispersion up to 10 times higher than predicted by Poisson models and highlight the impact of nonstationarity in assuming a constant mean rate within training-forecast intervals. The temporal rate variability is translated into a reduction of spatial precision by means of URZ models. First, counting processes are devised to capture rate distributions, considering the rate as a random variable. Second, we devise a data-driven method based on geodetic strain rate to spatially delimit the precision of URZs, assuming that strain/stress rate is related to the timescales of earthquake interactions. Finally, rate distributions are inferred from the available data within each URZ. We provide forecasts for the New Zealand National Seismic Hazard Model update, which can exhibit rates up to ten times higher in low-seismicity regions compared with SSMs. This study highlights the need to consider nonstationarity in seismicity models and underscores the importance of appropriate statistical descriptions of rate variability in probabilistic seismic hazard analysis.